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~qq~- Abstract. The following work presents a generalized (extended) finite element formulation for the 

advection-diffusion equation. Using enrichment functions that represent the exponential nature of 
CNJ ■ the exact solution, smooth numerical solutions are obtained for problems with steep gradients 

and high Peclet numbers (up to Pe = 25) in one and two-dimensions. As opposed to traditional 
stabilized methods that require the construction of stability parameters and stabilization terms, 



the present work avoids numerical instabilities by improving the classical Galerkin solution with 
an enrichment function. To contextualize this method among other stabilized methods, we show 
■ by decomposition of the solution (in a multiscale manner) an equivalence to both Galerkin/least- 

squares type methods and those that use bubble functions. This work also presents a strategy 
i for constructing the enrichment function for problems with complex geometries by employing a 

global- local approach. 
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The advection-diffusion equation represents transport by means of transmission and dissemina- 
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tion. The problem is significant both in its application to physical phenomena and as a precursor 
to studying more complicated equations like the Navier-Stokes equations. It is well known that the 
classical Galerkin method performs poorly for advection dominated transport problems. Spurious 
oscillations manifest in the solution due to the inherent truncation error of the discretized Galerkin 
approximation. Numerous strategies have been proposed in the literature to overcome this defi- 
ciency, which include: those that add artificial diffusion or employ non-centered discretization of 
the advection operator (also called upwind schemes) [U [21 El [H [5] , and multiscale models using 
bubble functions or wavelets [BJ [7]. In many cases, an equivalence can be constructed between 
these methods [3 [9]. For a complete review of stabilized methods for advection-diffusion, see [10] 
and the references therein. 
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Franca et al. [10] show that for the advection-diffusion problem with small k, small variations in 
the data produce large variations in the solution which represents the fundamental characteristic 
of instability. The discrete classical Galerkin formulation also inherits the lack of stability of the 
continuous problem. As the non-symmetric advection operator begins to dominate, the approxima- 
bility of the Galerkin finite element method degrades. For diffusion dominated flows, the solution is 
smooth, without sharp gradients or corners. In such cases, the underlying finite element subspace 
adequately captures the exact solution. Stabilized finite element methods improve the Galerkin 
formulation by reducing the truncation error (for example, Streamline-Upwind/Petrov-Galerkin) 
or increasing its approximability (for example, variational multiscale methods). 

Generalized (or extended) finite element methods (GFEM or xFEM) have been successfully 
applied in a number of contexts, including: solid mechanics (cracks) \12\ 113] . fluid-structure 
interaction (moving interfaces) [H], immiscible multiphase flows [15], and thermal problems with 
moving heat sources [16;. Other mesh free methods include the method of finite spheres [T7] and 
the reproducing kernel method [18]. This work presents the first application of the generalized 
finite element method to develop a stabilized- type finite element formulation. To formulate a 
robust stabilized formulation requires the construction of stabilization parameters and stabilization 
terms. In the generalized finite element context, one can improve the performance of the classical 
Galerkin formulation by designing an appropriate function that captures the solution better than 
the underlying finite element shape functions. 

We show in a subsequent section that there are many similarities between generalized and sta- 
bilized finite element methods. One can easily cast the generalized finite element method in such 
a way to obtain a stabilized formulation. Additionally, the structure of the resulting stability pa- 
rameter is identical to that of the Hughes variational multiscale method [6] . The methods differ in 
that, whereas stabilized methods aim to add additional terms to the formulation, the underlying 
Galerkin formulation need not change for the generalized finite element method. The methods also 
differ in that the generalized finite element method only enriches nodes in the vicinity of local phe- 
nomena such as a boundary layer. Stabilized methods are typically applied to the entire domain. 
Also, unlike stabilized methods, for generalized finite element methods the enforcement of Dirichlet 
boundary conditions is a bit more complicated as many of the enriched nodes fall on the boundary 
of the domain [19]. Another complexity of generalized finite element methods regards numerical 
integration techniques for mesh free methods |20j . 
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The primary aim of this work is to produce smooth solutions to the advection-diffusion equation 
by means of the generalized finite element method rather than a stabilized formulation. The nu- 
merical results presented show the successful application of the method to one and two-dimensional 
problems with high Peclet number (up to Pe = 25), with sharp corners, and boundary layer phe- 
nomena present. 

In the sections that follow, we first introduce the governing equations and the classical Galerkin 
formulation. We then discuss how the classical Galerkin formulation is enriched in the generalized 
finite element method. Since the generalized finite element method is heavily dependent on the 
enrichment function, we then present several enrichment functions appropriate for the advection- 
diffusion equation. Next, we discuss the enforcement of Dirichlet boundary conditions and the 
relationship between the generalized finite element method and stabilized methods. Lastly, we 
present some numerical examples in one and two-dimensions. 

2. GOVERNING EQUATIONS FOR ADVECTION-DIFFUSION 

Let be a bounded open domain, and r (T := Cl — f2, where f2 is the closure of Q) be its piecewise 
smooth boundary. Let a scalar unknown field be denoted by u : — ► M. Let the advection velocity 
field be denoted by a : — ► M. , where nd is the number of spatial dimensions and the diffusivity 
tensor be denoted as n : Q — ► ^ ndxnd . As usual, T is divided into two parts, denoted by T u and T*, 
such that r n n r* = and T n U T* = T. T u is the part of the boundary on which Dirichlet boundary 
conditions are prescribed, and T* is the part of the boundary on which Neumann (or diffusive flux) 
boundary conditions are prescribed. The governing equations for the advection-diffusion problem 
can be written as 

(1) ct • Vu - V • kVu = / in (! 

(2) u = vP on T u 

(3) kVu -n = t p on T* 

where vP is the prescribed Dirichlet boundary condition, t p is the prescribed flux, V is the gradient 
operator, / is the volumetric source term, and n is the unit outward normal vector to T. The fluid 
is assumed incompressible, (i.e. V • a = 0). 

The relative influence of advection and diffusion in equation ([1]) is expressed by the nondimen- 
sional global Peclet number, Pe = UL/k where L is the reference length, and U is a reference 
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velocity. In the finite element context, we define the element Peclet number as Pe h =|| a \\ Ii/2k, 
where || • || represents the L2 norm and h is the characteristic element length. 

2.1. Classical Galerkin formulation. We now define function spaces that will be used in the 
remainder of the paper. The functions spaces for the unknown scalar u(x) and its associated 
weighting function w(x) are 

(4) U :={u\u£ H 1 (£l),u = vP on T u } 

(5) W:={w\weH 1 (Q),w = Oonr u } 

The classical Galerkin formulation for the advection-diffusion problem is written as: Find u(x) £ U 
such that 

(6) a(w;u) + b(w;u) = l(w) V w £ W 
Let us define the bilinear forms as: 

(7) a(w;u) := / w(a ■ Vu) dSl 

Jn 

(8) b(w;u) := / Vw ■ kVu dO 

Jn 

and the linear functional as 

(9) l(w) := I wf dQ+ I wtdT 

Jn 7r* 

Once the weak formulation of the governing equations is established, the approximate solution 
based on the finite element method is determined in the usual manner. First one chooses the 
approximating finite element spaces, which (for a conforming formulation) will be the finite element 
function spaces for the unknown scalar u and the weighting function associated with u denoted by 
XA h C IA and W h C W respectively. The finite element formulation of the weak form reads: Find 
u h (x) £ U h such that 

(10) a(w h ;u h ) +b(w h ;u h ) =l(w h ) Vto'eW' 1 
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2.2. Stabilized methods. It is well known that the standard Galerkin method performs poorly 
for advection-dominated flow (i.e. Pe h S> 1) [21j . Excessive mesh refinement is necessary to avoid 
spurious oscillations that propagate throughout the solution. One popular method for dealing with 
such oscillations involves the use of stabilized methods (For examples of stabilized finite element 
methods for the advection-diffusion equation see |21| l22"j [23"1 [U [U EH [9] and the references therein) . 
Stabilized methods perturb the Galerkin formulation with mesh dependent terms. In some cases, 
the perturbation is designed in an ad hoc manner, while for other methods the perturbation pre- 
serves consistency by being residual based. In the present context, we introduce a stabilized- type 
formulation that is based on improving the classical Galerkin approximation by enrichment with a 
function that captures certain features of the solution. 

3. ENRICHED APPROXIMATION 

In this section we present a generalized finite element approach to solving the advection-diffusion 
equation by employing the partition of unity framework as detailed in [24l H2J ESI [26] . Let us divide 
the domain £l h into "Nele" non-overlapping subdomains fig (which in the finite element context 
will be elements) such that 

Nele 

(11) n h =\J n h e 

e=l 

The boundary of element Og is denoted by Tg . We then establish the standard finite element basis, 
Uq as the space of complete polynomials P fc (Og) of order < k over each element 

(12) := span{Ni}f =l where {N G 

where 2Vj is often referred to as the nodal shape functions, and F is the number of nodal finite 
element shape functions. In the numerical examples, piecewise linear approximations are used for 
N. Consider a set of overlapping arbitrarily shaped domains U{ centered at node i that define the 
support for each shape function, Ni. Also consider enrichment functions 7i m that represent special 
approximating properties in the proximity of local phenomena. A generalized finite element basis 
for U% that preserves a conforming approximation can be constructed as 

(13) U% := span {{N} U {AW m }£=i}ii 
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where E represents the number of enrichment functions for node i. The unknown scalar field u(x) 
may then be approximated as 



(a:) = ^ UiNi(x) + ^ N t (x) I ^«) m H m (i) j 

i i \ m / 



(14) u( 
where Ui represents the extra coefficients on the enriched nodes. 

Remark 1. Typically, when recovering the solution for a standard finite element analysis, the 
coefficients Ui are equal to the value of the unknown scalar field u{x) at the nodes. In the case of 
generalized finite elements, the coefficients Ui are not equal to the unknown scalar field at the nodes. 

Since our focus is on capturing the phenomena local to the boundary of the domain of interest, 
we need not enrich all the nodes in the domain. For the case in which the finite element space is 
enriched by only a single function 7i{x) at nodes in the local support of the boundary, the unknown 
scalar field u(x) may be approximated by 

(15) u(x) = ^2u i N i (x) + ^2u j N j (x)H(x) 

where I represents all nodes in the mesh and J represents the nodes that form a partition of unity 
for the function 7i{x) 

(16) J := {j € I : n rj + 0} 

To simplify the proceeding steps, we introduce the vectors of unknowns at each node, u and element 
shape functions N, 

(17) u T = [u u u' n ] ; N(x) = [JVi, • • • , N n , HN lt • • • , HN n ] 
where n is the number of nodes per element, such that 

(18) u{x) = N(x)u 

Note that for elements that do not have any enriched nodes, N = [Ni, ••• , -/V n ] and u T = 
[u~i, ■ ■ ■ , u n ] . After substitution of equation (|18p into equation (|10|) and noting the arbitrariness of 
w, the discrete form of the generalized finite element method of the advection-diffusion equations 
may be written as the assembled sum contribution of each element, 

Nele Nele 

(19) Aj2a(N;N)u + b(N;N)u = f(N) 

e=l e=l 



where A is the standard assembly operator. Equation (|19p is written in the familiar matrix form 

as 



(20) Ku = f 
where 

Nele 

(21) K = Aj2a(N;N)+b(N;N) 

e=l 
Nele 

(22) / = A £ /(AT) 

e=l 

Equation (|20p represents a linear system of equations with the number of unknowns equal to the 
number of nodes plus the number of enriched nodes. 

4. DESIGN OF ENRICHMENT FOR ADVECTION-DIFFUSION 

Typically, the enrichment function or set of functions 7i{x) are chosen to provide an improved 
estimate for the solution in the vicinity of local phenomena. In some cases, the exact solution is 
known in the immediate vicinity or at least the nature of the exact solution is known. In such cases, 
Ti(x) may be choses as a particular part of the solution. For complex geometries, or for problems 
for which little is known about the solution, we introduce a global-local approach to choosing the 
enrichment function. In a subsequent section, we explore the dependence of the solution on the 
enrichment function selected. 

4.1. One-dimensional problem. Consider the advection-diffusion problem in one dimension on 
a domain of unit length with the following boundary conditions 

, . du d 2 u ^ 

(23) a- k— 7t = 1 in O 

ox ox 1 

(24) u(x) = on r 

Assuming a constant advection speed, a, and diffusivity, n > 0, the exact solution to this problem 
is 

1 ( 1 — e ax l K \ 

(25) u(x) = -\x 



a \ 1 - e a / K 
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The nature of the solution is highly dependent on the Peclet number, Pe h . For high Pe h a thin 
boundary layer is present toward the outflow boundary. The exact solution for increasing Pe h is 
shown in Figure [TJ The simplest enrichment function, motivated by the exact solution, is to use 

(26) H a (x) = e>* 

where 7 = ol/k. In this case, we take advantage of the exponential nature of the solution. Anal- 
ogously, in linear fracture mechanics problems, one knows that the stress field grows as 1/y/r, 
where r is the radius from the crack tip. For advection-diffusion we know that the solution has an 
exponential nature. 

Of primary interest for flow problems is the application of the generalized finite method to cap- 
turing the localized nature in the boundary layer. As such, the enforcement of Dirichlet boundary 
conditions becomes more involved as many of the enriched nodes lie on the Dirichlet boundary. 
The enforcement of Dirichlet boundary conditions is discussed in detail in the next section, but we 
wish to construct here an alternative enrichment function that does not interfere with the Dirichlet 
boundary. Building on the enrichment function 7i a {x) as motivated by the exact solution, we can 
construct an alternate normalized enrichment function 7ib(x) that vanishes on the boundary, 

(27) H b (x) = 1 - 

where again 7 = ol/k. Figure [2] shows TCb(x) for various Pe. 

Closer examination of equation (|27p reveals that 7ib{x) is not well-defined numerically for 7 = 
or 7 much greater than 800, which represent the cases of pure diffusion or advection dominated flows 
respectively. An alternative enrichment TC c (x) may be selected that is a polynomial representation 
of the exponential nature of TCb(x) 

(28) H c {x) = 1 - x 7 

TC c (x) both vanishes at the outflow boundary and is well-defined for pure diffusion and advection 
dominated flows. 

4.2. Global-local enrichment function. In multi-dimensional problems on complex geometries, 
choosing an appropriate enrichment function becomes more difficult. As such, one can construct 
an enrichment function based on the solution to the problem at low Pe h , for which the classical 
Galerkin formulation is stable. In the global-local approach, one begins with a stable solution at 
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low Pe h and iteratively increases the Pe h , using the previous solution as the enrichment function. 
The algorithm goes as follows, 

Global-local solution strategy 

Solve the Galerkin problem (equation (jlOp ) for Uq, with Pe^ = 1.0 
n = 4 

APe h = (Pe 1 ^ — Pefy/N, where N is the number of iterations 
for i = 1 to iV 



The preceding solution strategy may also be considered as a continuation method in the con- 
text of flow problems [27]. The aim of a continuation method is to solve the problem for a high 
characteristic flow speed using information from the solution at a lower characteristic speed. 



Among the various methods available in mesh-free methods for the enforcement of Dirichlet 
boundary conditions, one of the simplest ways to implement is the penalty method. The penalty 
method introduces no additional unknowns and preserves the banded structure, positive definite- 
ness, and symmetry of the stiffness matrix. 

5.1. Penalty enforcement. We begin with a variational statement from which the advection- 
diffusion equation may easily be derived in one dimension. Consider the variational functional 



where V is the gradient operator, || • || represents the appropriate norm, and / represents the 
body force. Solving the advection-diffusion equation is equivalent to extremizing H(u) subject to 
u(x) = uq on T M . Applying the stationarity condition to H(u), 



Hi = nti, Pef = Pe^i + APe h 

Solve the GFEM problem (equation flU) for u 4 h , with Pe = Pef 



end 



5. ENFORCEMENT OF DIRICHLET BOUNDARY CONDITIONS 



n(u) 



(29) 




(30) 
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we have 

f du d^u 

(31) / a— - k— — jr - / dO = ; s.t. u(x) = n on T u 

Jn Ox ox 

Using a penalty method, we change the constrained extremization problem to an unconstrained 
extremization problem by adding a penalty functional, AA(n)/2 to n(n) 

(32) A(n) = / ~{u{x) -n ) 2 dr 



where A is the penalty parameter. Extermizing the sum of both functionals n(n) and A(n), the 
solution will satisfy the Dirichlet boundary conditions as A approaches oo. The resulting functional, 
that incorporates a weak enforcement of the essential boundary conditions is written: 

(33) / a— - k— -fdQ + X / (u(x) - n ) dr = 

Jn ox ox J r u 

For mesh-free methods, the enforcement of essential boundary conditions must be handled with 
great care to avoid over-constraint phenomena. To eliminate spurious results when enforcing the 
Dirichlet boundary conditions weakly, one can use trapezoidal numerical integration rather than 
Gauss integration to evaluate the boundary terms (/ r „(-) dr) in equation (|33|) [19] . When en- 
richments are used that do not vanish on the Dirichlet boundary (for example TC a ), the boundary 
conditions will be enforced in a weak fashion according to equation (|33|) . 



6. GFEM AND STABILIZED METHODS 

Due to the proliferation of stabilized methods for the advection-diffusion problem, as pointed 
out in the introduction, we wish to cast the generalized finite element formulation for advection- 
diffusion in a different manner to explore its relationship with stabilized methods. To explore this 
relationship, we begin by establishing a variational mulitiscale framework, in which we separate 
the problem into two sub-problems: the first representing the standard finite element solution, and 
the second the enriched finite element solution. We then solve the enriched problem in terms of 
the standard variables and substitute back into the standard finite element problem. The result is 
a stabilized formulation equivalent to Galer kin/least-squares (GLS) or streamline- upwind/Petrov- 
Galerkin (SUPG) and the structure of the stabilization parameter is identical to that of the Hughes 
variational multiscale method. An equivalence between stabilized methods such as GLS and multi- 
scale methods that employ bubble functions has previously been developed in works such as [HI [9] , 
but we wish to extend this equivalence to the generalized finite element context. 
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Consider the advection-diffusion equation enriched with a single function, 7i{x) and constant 
diffusivity k. We first decompose the unknown field in equation (|15p and its weighting function 
into a standard Galerkin part, u, and an enriched part, u' , as follows: 

(34) u{x) = u{x) + u (x) ; w(x) = w(x) + w'(x) 
where 

(35) u(x) = ^2uiNi(x) = Nu ; u{x) = ^u' j N j (x)H{x) = N'u 

iei jeJ 

(36) w(x) = ^2wiNi(x) = Nw ; w'(x) = ^w' j N j {x)H{x) = N'w' 

iei jeJ 

The weak form (equation (jlOp ) may now be written as 

(37) a(w + w';u + v!) + b(w + w';u + v!) = l(w + w') 

Due to the linearity of the weighting function and the trial solution, we can decompose (|37p into 
two subproblems: the Galerkin subproblem, 

(38) a(w; u) + a(w; v!) + b(w; u) + b(w; u') = l(w) 
and the enriched subproblem 

(39) a(w'; u) + a(w'; u') + b(w'; u) + b(w ; u) = l(w') 
Using equation §5§ and integration by parts as follows, 

(40) / w(a ■ Vu) dtt = - / a ■ Vw ■ u dQ = c(w; u) 

Jn Jn 

(41) / !«V • kVu dO = - / wkV 2 u dn = d(w; u) 

Jn Jn 

The Galerkin subproblem and enriched subproblem may alternatively be written as 

(42) a(w; u) + b(w; u) + c(w; u') + d(w; u') = l(w) 
and 

(43) a(w'; u) + a(w'; u) + b(w'; u) + d(w'; u) = l(w') 
respectively. 

li 



Consider the enriched subproblem, equation (|43p . which can be written as 

(44) / w'(a ■ Vu) + Vw' ■ kVu = - / w'r dO 

Jn Jn 

where f = —ot-Vu + kV 2 u- / '. Using equations ([55]) and ([55]) . the above equation can be expressed 

in discrete form as 

(45) w' T [ N'a T X7N' + kVN' t VN' dfiu' = -w' T [ N ,T r dO 

in Jn 

Noting the arbitrariness of w', we have 

(AR) ii = J J± 

1 ' f Q N' T a T \7 N' + K VN' T VN' dQ 

Due to equation (|35p the enriched solution can be written as a function of the variables of the 
Galerkin subproblem as 

N ' lo N ' T dfif 

47 U ( X ) = —r W r —> 7r / = ~ Tf 

f n N a T V N + kV N V N dQ 
where we have made the assumption that f is constant over the element domain, which is exactly 
true for a constant body force using linear elements, and we have introduced the familiar intrinsic 
time scale, r, also known as the stabilization parameter 

N'f n N' T dft 

(48) T ~ / n iV /T a T VAr' + K VN' T VN' dn 

Substituting equation ([15]) into the Galerkin subproblem, equation (|4*2|) . we have 



(49) a(w; u) + b(w; u) + c(w; rr) + d(w; rr) = l(w) 

Using the general inner product, (w;u) = f n (-)d£l Equation (|4^|) can be written more explicitly as 

(50) {w, a ■ Vm) + (Vw, kSJu) + (a-Vw + kV 2 w, r(a • Vn - kV 2 m)) = 

(51) (w,f) + (T(a-Vw + KV 2 w),f) 

which is precisely the stabilized formulation proposed in [28] . In [28j , the stabilization parameter is 
defined based on error analysis considerations. In the present work, we have derived a stabilization 
parameter based on a decomposition of the generalized finite element method. To further illustrate 
the relationship between the generalized finite element method of advection-diffusion, consider the 
structure of r in equation (|48p . which is identical to the structure of the stabilization parameter 



derived using the Hughes variational multiscale method in [B]. In [BJ, for N' , the authors employ 
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bubble functions that vanish on the element boundary. The nature of the bubble functions preserves 
convergence as the bubble functions do not interfere with the nodal values. In the generalized finite 
element framework, N' does not vanish on the element boundary, but is defined over the entire 
domain. In this case, convergence is ensured based on the qualities of the partition of unity [24] . 

Clearly, casting the generalized finite element method into a variational multiscale framework 
elucidates its apparent similarities with stabilized methods. 



In this section, we present some numerical examples that illustrate the potential of the generalized 
finite element method for producing stable solutions to the advection-diffusion equation for high 
Pe h . Due to the sharp corners local to the enriched elements, a 100-point gauss integration scheme 
is used to integrate the stiffness terms in equation (|2ip for any element with enriched nodes. Also, 
for the examples that enforce the boundary conditions weakly, the boundary integrals are integrated 
using two-point trapezoidal numerical integration. 

7.1. Linear elements in one-dimension. The first example we present is the advection-diffusion 
equation in one-dimension, with a constant advection speed a and diffusivity k, and a unit source. 



The problem is solved on a unit domain discretized into six equal length elements. The results are 
shown in Figure Notice that the generalized finite element method performs well for high Pe. 
Also notice that the solution between the nodes is indistinguishable from the exact solution for the 
generalized finite element method. The results shown in Figure were obtained by enriching the 
last two nodes with enrichment function Tib • To show the equivalence of enriching with Ti a , Tib, or 
Ti c a comparison is shown in Figure HI Since enrichment functions Ti a and Ti c do not automatically 
vanish on the boundary, the boundary conditions are enforced weakly in those cases. 

7.2. Linear elements in two-dimensions. The same enrichment function used in one-dimension 
is also applicable in two-dimensions. To illustrate this point, we present the results from the 
advection-diffusion equation applied to a unit square domain, with constant advection speed and 
diffusivity, and a unit source term. The mesh used for this example is shown in Figure [5] along with 



7. NUMERICAL RESULTS 



(52) 



du d 2 u 




(53) 



u = on r 
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the location of the enriched nodes. The results shown in Figure [6] (b) were produced by enriching 
with the following function 

1- eY x i- e iv 

(54) = (i-— rXi-j— T ) 

with 7 = a/ k. Figure [6] (a) shows the standard Galerkin solution to the problem, which is highly 
oscillatory. Figure [U] (c) shows a comparison with results obtained by using the standard Galerkin 
method with exponential shape functions rather than linear elements. 

7.3. Linear elements in one and two-dimensions using the global-local approach. Be- 
cause the proper enrichment function is not always an obvious choice for every problem, we present 
some numerical examples for which the solution to the problem at low Pe h is used for the enrich- 
ment function. Figure [7] shows the results from using the global-local approach to solve the same 
ID problem described above with a unit source, vanishing u on the boundary, on a unit domain. 
Whereas using enrichment function 7i.b produces accurate results for high Pe h (as shown in Figure 
E]), the global- local approach provides accurate solutions to Pe h ~ 3. In two-dimensions, a smooth 
solution is obtained for Pe h = 7, shown in Figure [HJ 



7.4. Thermal boundary layer. The thermal boundary layer problem was first analyzed in 
and represents the simulation of fully developed flow between two plates, the top plate moving with 
a unit velocity in the x-direction. The domain and boundary conditions are shown in Figure [9j The 
highest global Pe for this problem is ~ 700, whereas the element Pe h using the elements nearest 
the top plate is 25. The mesh with enriched nodes along the outflow boundary is shown in Figure 
[TOl As shown in Figure [TTT smooth solutions are obtained using enrichment function TC C , enforcing 
the boundary conditions weakly. Strongly enforcing the boundary conditions causes spurious oscil- 
lations as predicted in [19]. Notice that the Galerkin solution has wiggles that propagate through 
the entire domain. The results compare well with stabilized methods presented in |22[ [28] for the 
SUPG method, and [6] for GLS and the Hughes variational multiscale method, further illustrating 
the close relationship between stabilized and generalized finite element methods. 

8. CONCLUSIONS 

We have introduced a generalized finite element formulation for the advection-diffusion equation 
with stabilizing properties. Since the methods is highly dependent on the enrichment function 
selected, we have provided various enrichment functions with different qualities and shown an 
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equivalence between them. We have also explored the relationship between generalized finite ele- 
ment methods and stabilized methods. The similarity between GLS-type methods and stabilized 
methods using bubble functions has been extended to the generalized finite element context by 
decomposing the problem in a variational multiscale manner. The numerical examples presented 
show that the method is applicable to high Pe flows and works for problems with non-constant 
boundary conditions. The examples also illustrate the method's similarity to stabilized methods. 
Finally, this work emphasizes the need for a broader interpretation of stabilized methods for flow 
problems including a more comprehensive theory from which stabilized and generalized methods 
emanate from. 
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Figure 1. One-dimensional advection-diffusion: exact solution for various Pe. 
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Figure 2. ID enrichment function as a function of the radius r = 1 — x from the 
outflow boundary. 
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3. Finite element solutions to the advection diffusion equation for various 
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FIGURE 4. Generalized finite element solutions to the advection diffusion equation 
for Pe h = 5. 
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Figure 5. Finite element mesh for the 2D constant coefficients problem showing 
the enriched, partially enriched, and unenriched elements. 
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(c) 

Figure 6. Results for 2D advection diffusion equation for Pe h = 35 (a) Galerkin 
finite element solution (b) generalized finite element solution using Tt^ (b) Galerkin 
solution using exponential finite elements. 
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Figure 7. Generalized finite element solution to the ID advection-diffusion equa- 
tion using global-local enrichment for various Pe h . 
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Figure 8 . 2D constant coefficients on a regular domain: (a) Galerkin finite element 
solution (b) generalized finite element solution using global-local approach, enforcing 
Dirichlet boundary conditions weakly, Pe h = 7, A = 10 6 . 
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Figure 9. Thermal boundary layer problem: geometry and boundary conditions. 
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Figure 10. Finite element mesh for the thermal boundary layer problem showing 
the enriched, partially enriched, and unenriched elements. 
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(a) 



(b) 




(c) (d) 

Figure 11. Thermal boundary layer: (a) Galerkin finite element solution Pe h = 
12.5 (b) generalized finite element solution using 7i c , enforcing Dirichlet bound- 
ary conditions weakly with A = 10 10 , Pe h = 2.5 (c) Galerkin finite element solution 
Pe h = 25 (d) generalized finite element solution using 7i c , enforcing Dirichlet bound- 
ary conditions weakly with A = 10 10 , Pe h = 25. 



21 



